from pulp import *

NODE_COUNT = 28
SOURCE_NODE = 27
SINK_NODE = 19

#capacity and price
cs = [[0 for x in xrange(NODE_COUNT)] for x in xrange(NODE_COUNT)]
a = [[0 for x in xrange(NODE_COUNT)] for x in xrange(NODE_COUNT)]

def connect(i, j, capacity, c):
    cs[i][j] = cs[j][i] = capacity
    # add the cost
    other = range(20,27)
    if (i in other and j not in other) or (j in other and i not in other):
        a[j][i] = 1
        a[i][j] = 1

#define edges
connect(0	,	7	,	10	, cs)
connect(0	,	6	,	20	, cs)
connect(1	,	7	,	10	, cs)
connect(1	,	8	,	10	, cs)
connect(2	,	3	,	15	, cs)
connect(2	,	8	,	40	, cs)
connect(2	,	9	,	25	, cs)
connect(3	,	9	,	5	, cs)
connect(4	,	9	,	10	, cs)
connect(4	,	10	,	10	, cs)
connect(5	,	10	,	10	, cs)
connect(6	,	7	,	10	, cs)
connect(6	,	18	,	20	, cs)
connect(7	,	8	,	10	, cs)
connect(7	,	18	,	10	, cs)
connect(8	,	12	,	10	, cs)
connect(9	,	11	,	15	, cs)
connect(9	,	12	,	20	, cs)
connect(10	,	11	,	20	, cs)
connect(11	,	13	,	20	, cs)
connect(11	,	14	,	20	, cs)
connect(12	,	15	,	30	, cs)
connect(12	,	17	,	30	, cs)
connect(12	,	23	,	10	, cs)
connect(13	,	14	,	10	, cs)
connect(13	,	15	,	40	, cs)
connect(14	,	15	,	10	, cs)
connect(14	,	19	,	30	, cs)
connect(15	,	16	,	10	, cs)
connect(15	,	24	,	10	, cs)
connect(15	,	19	,	30	, cs)
connect(16	,	19	,	50	, cs)
connect(16	,	24	,	200	, cs)
connect(16	,	25	,	200	, cs)
connect(17	,	18	,	30	, cs)
connect(17	,	23	,	10	, cs)
connect(18	,	20	,	10	, cs)
connect(18	,	21	,	10	, cs)
connect(20	,	22	,	50	, cs)
connect(20	,	26	,	50	, cs)
connect(21	,	22	,	50	, cs)
connect(21	,	23	,	50	, cs)
connect(22	,	23	,	50	, cs)
connect(22	,	24	,	50	, cs)
connect(22	,	25	,	50	, cs)
connect(22	,	26	,	50	, cs)
connect(23	,	24	,	50	, cs)
connect(25	,	26	,	50	, cs)
#super source edges:
connect(0, SOURCE_NODE, sys.maxint, cs)
connect(1, SOURCE_NODE, sys.maxint, cs)
connect(2, SOURCE_NODE, sys.maxint, cs)
connect(3, SOURCE_NODE, sys.maxint, cs)
connect(4, SOURCE_NODE, sys.maxint, cs)
connect(5, SOURCE_NODE, sys.maxint, cs)

# create the max flow problem
maxFlow = LpProblem("Max Flow Problem",LpMaximize)
# create the lowest price problem
minCost = LpProblem("Min-Cost Flow Problem",LpMinimize)

# create variables with constraints
f = [[] for i in range(NODE_COUNT)]
for u in range(NODE_COUNT):
    for v in range(NODE_COUNT):
        f[u].append(LpVariable("flow_"+str(u)+"_"+str(v),0,cs[u][v],LpInteger))

# add the objective
maxFlow += lpSum([f[SOURCE_NODE][i] - f[i][SOURCE_NODE] for i in range(NODE_COUNT)]), "Max Flow"
ob = []
for u in range(NODE_COUNT):
    for v in range(NODE_COUNT):
        if a[u][v] != 0:
            ob.append(a[u][v]*f[u][v])
minCost += lpSum(ob), "Min-Cost Flow"

# add flow conservation
for u in range(NODE_COUNT):
    if u != SOURCE_NODE and u != SINK_NODE:
        maxFlow += lpSum([f[v][u] for v in range(NODE_COUNT) if v != u]) <= lpSum([f[u][v] for v in range(NODE_COUNT) if v != u]), "(1)flow conservation for %s" % u
        maxFlow += lpSum([f[v][u] for v in range(NODE_COUNT) if v != u]) >= lpSum([f[u][v] for v in range(NODE_COUNT) if v != u]), "(2)flow conservation for %s" % u
        minCost += lpSum([f[v][u] for v in range(NODE_COUNT) if v != u]) <= lpSum([f[u][v] for v in range(NODE_COUNT) if v != u]), "(1)flow conservation for %s" % u
        minCost += lpSum([f[v][u] for v in range(NODE_COUNT) if v != u]) >= lpSum([f[u][v] for v in range(NODE_COUNT) if v != u]), "(2)flow conservation for %s" % u

# solve maxflow
maxFlow.solve()

print "Status:", LpStatus[maxFlow.status]

##for variable in maxFlow.variables():
##    print variable.name, "=", variable.varValue

print "Maximum flow = ", value(maxFlow.objective)

# add the flow constraint for minCost
minCost += lpSum([f[SOURCE_NODE][i] - f[i][SOURCE_NODE] for i in range(NODE_COUNT)]) <= value(maxFlow.objective), "(1)Flow Constraint"
minCost += lpSum([f[SOURCE_NODE][i] - f[i][SOURCE_NODE] for i in range(NODE_COUNT)]) >= value(maxFlow.objective), "(2)Flow Constraint"

# solve mincost
minCost.solve()

print "Min Cost Status:", LpStatus[minCost.status]
print "Min-Cost flow = ", value(minCost.objective)